Perturbative expansion using variational phonon basis for Holstein model 
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A simple variational displacement phonon basis, obtained through the modified Lang-Firsov 
(MLF) transformation, is proposed to study the Holstein model. This phonon basis contains only 
one variational parameter, but capable of describing lattice distortions at distant sites from the 
charge carrier. Perturbation method based on this MLF basis is employed to calculate the single- 
electron ground-state energy and the dispersion of the polaronic band. The ground-state (fe=0) 
energy obtained up to the second-order perturbation within this approach agrees well with the 
available numerical results for the entire range of coupling strength. 
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Self-localization of electrons in crystal lattices due to 
interaction with phonons constitutes one of the central 
problems in condensed matter physics. The simplest Hol- 
stein model describing coupling of conduction elec- 
tron with dispersionlcss phonons, has been studied in 
weak H] as well as strong coupling [|| limits following ei- 
ther variational ^, |j or perturbation methods ||, ^ . 
Numerical approaches including exact diagonalization 
(ED) techniques [l^. Quantum Monte Carlo calcu- 
lations Global-Local (GL) method density- 
matrix-renormalization-group (DMRG) technique 
and recently developed ED technique with variational 
Hilbert space are quite successful in evaluating 
the ground-state properties including various correlation 
functions in the entire range of electron-phonon (e-ph) 
coupling. The analytical methods, however, suffer from 
many limitations. Migdal approximation |^ is valid for a 
weakly-coupled adiabatic system while Lang-Firsov (LF) 
canonical transformation H is suitable only in the strong- 
coupling limit. No common consensus regarding proper 
choice of the phonon basis and applicability of pertur- 
bation theory has yet been reached for the whole range 
of e-ph coupling. So far the perturbative method has 
been employed using the LF basis only and no perturba- 
tion calculation based on the modified LF (MLF) basis 
is available. Recently, we [D investigated the conver- 
gence of the perturbation expansion on a two-site sin- 
gle electron Holstein model using different phonon bases 
obtained through LF, MLF and MLF with squeezing 
(MLFS) transformations. It has been observed that for a 
wide region of e-ph coupling the perturbation corrections 
within the MLF method are much smaller and the con- 
vergence is much better compared to the LF approach. 

Inspired with this result, we now attempt to apply the 
perturbation method based on a properly chosen MLF 
basis to the one-dimensional (infinite site) Holstein model 
and evaluate the perturbation corrections to the ground- 
state energy, wave function and the dispersion. The Hol- 



stein Hamiltonian for a single electron reads as 

i i i,5\ 



ujg 



(1) 



where i =1 to N, denotes the site. Ci (c|) is the annihila- 
tion (creation) operator for the electron at site i and rii 
{—clci) is the corresponding number operator, g denotes 
the on-site e-ph coupling strength and t is the usual hop- 
ping integral, bi and 6| are the annihilation and creation 
operators, respectively, for the phonons corresponding to 
interatomic vibration at site i and uj is the phonon fre- 
quency. Si runs over nearest-neighbor sites. 

The spread and depth of the polaron may be studied 
conveniently within the framework of the MLF trans- 
formation where the lattice deformations produced by a 
charge carrier at different neighboring sites are treated as 
variational parameters ||^, ^. The MLF transformation 
on Eq. (1) gives H = e^He^^ 
where 



R = ^n,[Ao(6l-fe,) + Ai^(6t^,^-&,+5j 

i 8i 

+ A2 ^(&I+52 - + - + Ai ^ (&!+,-, - + 

&2 Si 

and Aq, Ai, A2, ... represent the variational lattice defor- 
mations around the charge carrier produced at the same, 
first-, second-, etc. neighbouring sites, respectively. 5i 
runs over the I th nearest neighbour sites. For LF trans- 
formation, Xo—g and Ai = A2 = ... = A; = ... = 0. 

The MLF transformed Hamiltonian is then obtained 
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H = uj^blbi + ep^Ui + Ht + 

i i 

where 

i,Si 

Hs = LLi^{bl +bi)[(g - Xo)n^- Xi^n^+Si 

i Si 

<53 

= Ao(6l-6.) + Ai^(feU-6.+5j 

<5i 



(2) 



and 



ep = e - (2.g - Aq) Aq + 2tj (Af + A^ + A^ + ...) 

is the polaron self-energy. 

For the perturbation calculation following Marsiglio 0] 
we choose the basis set for the one-electron case as 



cj\0)e\ni 712 .. nN)ph 



(3) 



which represents a state with Zth site occupied by the 
electron and ni,n2, ...,rtAr being the phonon numbers at 
sites 1, 2, ... , N, respectively in the MLF displaced 
phonon basis. I takes values from 1 to while ni,ri2, 
... from zero to infinity covering the whole Hilbert space. 
Ht and Hs are non-diagonal and the remaining part of 
the Hamiltonian is diagonal in this representation. The 
unperturbed energy of the state, described in Eq. (3), is 
given by -Sj^^. j = ep+nx <^ where nT{= '^0 represents 
the total number of phonons in that state. The ground 
state corresponds to tt-t = and has A^-fold degeneracy 
since the electron can occupy any one of the N sites. This 
degeneracy is lifted by considering the first-order energy 
correction which results in a polaronic band with energy 
Ek = £p — 2te cos(fca) where tg — <exp(— {(Aq — Ai)^ -I- 

(Al-A2)2 + (A2-A3)2 + ..}). 

The eigenstates are now given by 



-^^exp(ifc.^;) 



4|0)e|0) 



ph 



(4) 



and the unperturbed ground state \ipG^) corresponds to 
k = 0. 

For perturbation calculation it is very tedious and in- 
convenient to deal with a phonon basis containing too 
many variational parameters. Previously we observed 
that a choice of variational parameters which minimizes 
the ground-state energy, obtained within MLF and zero- 



phonon averaging (ZPA) (i.e., from ph{0\H\0)pfi) leads 
to a phonon basis where the convergence of the per- 
turbation series is satisfactory. Here minimization of 



Ek=o = — (equivalent to the ground-state energy 
obtained by ZPA of H) with respect to parameters (Aq, 
Ai, ...,Ai) gives the general recursion relation between the 
variational parameters 



A;-i A(_2 — 2A/_i + A/ 



(5) 



A/ A(_i — 2A; + A;+i 

which is satisfied by a simple choice for the parameters 



Ai _ A2 _ A3 
Ao Ai A2 



Xi-i 



(6) 



A sum rule is also followed by A;'s; for one-dimensional 
chain this is given by 



Ao + 2Ai + 2A2 + ... + 2A, + .. 



(7) 



Eq. (7) is a consequence of the fact that the in-phase 
phonon mode (zero momentum) being a sum of lattice 

displacement operators at different sites {J^ii^i + ^!)) 
couples only with the total number of electrons, which 
is a constant of motion. Eqs. (6) and (7) allow one to 
express all A; 's in terms of a single variational parameter 
r. Thus, the obtained MLF phonon basis, though being 
capable of considering lattice distortions at all sites of the 
infinite chain, contains only one variational parameter. 
The polaron self-energy may then be written in terms 

of r as Cj, = e — ojg^ ''^~''(i^^3~'~'' ^-nd the effective hop- 
ping te = t exp{—g^x^) where x = 

The expression for the off-diagonal matrix element be- 
tween the ground state and an excited state is given by 



H, 



'nn},G = {Mn^}\{Ht+HsM 



G I 



N 



^J{...ni-i\ nil ni+i\....) 



[(1 - x)6niS - X[r{6n,_i,l + Sni + ul) 



(8) 



where?/ = gx{l — r), 

UL = + n;_2 + ?J;-3 + 

riR = + n/+2 + n;+3 + 

A = n;+i + n/_2 + 2(n/+2 4- n/_3) 

+ 3(n;+3 + n/_4) + .., 

B = ni-i + ni+2 + 2(n/_2 + ni+3) 

+ 3(n,_3+n,+4) + .-. (9) 

The second-order correction to the ground-state energy 
is given by 
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— [-r-<,(l-r) ] - 



2 (3 + r) 



(l-hr)3 

1 ^ ^ t2y2«^ [(-1)"^ + (-1)"^ 
^ / n7>2 {.■■■ni-i\ ni\ ni+i\....) 



(10) 
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and the first-order correction to the ground-state wave 
function is obtained as, 

1 ^ 

m p 

N jj' 

Now one has to make a proper choice of r. In addition 
to Eq. (6) the minimization of Ek=a yields a condition 

Lor = tc{l — r)'^ , which gives 

r = 

It is interesting to note that for this choice of r the off- 
diagonal matrix element between the ground state and 
any excited state with a single phonon (jit = 1) becomes 
zero and this is consistent with our previous studies 
Correspondingly, the first term in Eq. (10) as well as in 
Eq. (11) vanishes. 

3 I 1 1 1 — , 1 1 1 1 1 
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FIG. 1: Plot of logarithm of the effective mass m* /mo as 
function of g for t = uj = 1.0 within MLF and LF approach 
using the unperturbed ground-state wave function. '4-' shows 
the results of Ref [14]. 



For real materials of interest the Hamiltonian may involve 
other interactions in addition to the Holstein-type e-ph 



interaction and the perturbation calculations for such 
complicated systems may not be possible. The common 
practice in such cases is to make LF or MLF transfor- 
mation followed by ZPA. It may be mentioned that ZPA 
is equivalent to the zeroth (first) order of perturbation 
as far as the wave function (energy) is concerned. The 
success of the ZPA depends on the choice of the phonon 
basis. One of the advantages of the MLF phonon ba- 
sis over LF basis is that the second-order perturbation 
correction (£^q ' ) within MLF basis is much smaller than 
that within LF basis except in the strong-coupling regime 
where both the bases become equivalent. For instance, 
for w = t 1, E^^\mLF) = -0.01746 and -0.1849 

while e''^\lF) = -1.3538 and -0.4801 for g = 1.0 and 
1.7, respectively. 

A comparison of the effective mass of the polaron ob- 
tained within the MLF and LF with ZPA may be help- 
ful. We have calculated the effective mass of the polaron 

using the standard formula ^ — ^ ^ fj^P \k=o (where 
mo = ^) within the ZPA of LF and MLF Hamiltonian. 
These results along with the exact numerical result of 
Ref. Ilj] are presented in Fig. 1. The proximity of the 
MLF results with the exact one indicates that the MLF 
basis is a better choice. 
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FIG. 2; Plot of the ground-state energy (in unit of t) ob- 
tained by considering perturbation corrections up to second- 
order within MLF, LF method and also within the weak- 
coupling (self-consistent Migdal approximation) and strong- 
coupling perturbation expansion (SCPT) obtained from Ref. 
[7] as a function of g for t = u) = 1.0. 
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In Fig. 2 we have shown the ground-state energy cal- 
culated up to the second-order perturbation correction 
within the LF (corresponds to r = 0) and the MLF 
phonon basis. Energies calculated from the weak cou- 
pling (self-consistent Migdal approximation) and strong 
coupling perturbation limit, obtained from (Eqs. (7) and 
(C.17) of) Ref. j7|, arc also plotted. From comparison 
with the GL results we find that the MLF energy al- 
most coincides with the exact energy in the entire regime 
of g. The weak-coupling method fails for large g whereas 
the strong-coupling perturbation expansion of Ref. |^ 
breaks down and the LF perturbation result deviates 
much from the exact energy for low g (< 1.4). 

To find out the energy dispersion within the MLF 
method the energy i?^. has been calculated up to the 
second-order perturbation. The second-order perturba- 
tion correction to within the MLF approach is ob- 
tained as 



UJ 



(3 + 
(l + r)3 



[ujr — te(l — r)^cos(fca) )^ 



-tg sm [ka) 



(12) 



(1 + r) 

1 „ „ t2y2nr (^2 + Y^ + 2XYcos{2ka)) 

N 2^ 2^ 



^ ^ ("^'^) (■■•"-;-!! nil ni+i\....) 



where X = (-l)"^r^ and Y = (-l)"«r^ 
A and B are defined in Eq. (9). 

The corresponding expression for i?^' 
approach is 



and 



(2) 



within the LF 



2il 

UJ 



°° 2 0^" 

rcos^(fca) 

n=l 



oo oo 

+ EE 



2(n+m) 



(n -t- m) n\ m\ 



(13) 



where tg = t e ^ 



(2) 

We have calculated hence Ek up to the second- 

order perturbation within the LF and MLF approach. 

(2) (2) 

The second-order correction Ej. and the ratio E^. /Ek 
for different k values are listed in Table I for several values 
of g for both the MLF and the LF approach for t = w = 1. 

The following features are evident from our study: 
(i) The second-order MLF corrections are appreciably 
smaller than the corresponding LF corrections for low 
values of g and k. (ii) For small and intermediate values 
of g (e.g. g — 0.5 and 1.0) the MLF correction increases 
with increasing fc, however the MLF corrections are less 
than the corresponding LF values in the entire region of 
k except at fca = 7r/2 where both the MLF and LF cor- 
rections are same. In fact, for ka = tt/2 the value of r 
becomes zero within the present MLF approach so the 
MLF basis reduces to the LF basis for ka = tt/2. 



(2) 

TABLE I: ka/iT, second-order correction energy {E). , in unit 



of t) and realtive second-order correction {e'/^^ / Ek) for MLF 
and LF method respectively, (for g = 0.5, 1 and 2.0, t = ui — 
1.0) for the ground state of many-site one polaron problem. 



MLF 



LF 



ka/-K 



E 



(2) 



Ei 



0.5 



0.00 
0.10 
0.20 
0.30 
0.40 
0.50 



-0.0010180 0.000482 -1.014952 0.359587 

-0.0094835 0.004680 -0.953203 0.355067 

-0.0364071 0.020493 -0.791541 0.343899 

-0.0881567 0.062911 -0.591716 0.336728 

-0.1841571 0.189624 -0.430054 0.370296 

-0.3683048 0.595669 -0.368305 0.595669 



1.0 
0.00 
0.10 
0.20 
0.30 
0.40 
0.50 



-0.0174239 
-0.0501885 
-0.1526451 
-0.3373983 
-0.5934159 
-0.6403983 



0.007045 
0.020757 
0.067293 
0.162439 
0.310214 
0.390392 



-1.353833 
-1.285706 
-1.107347 
-0.886884 
-0.708525 
-0.640398 



0.438192 
0.430657 
0.409736 
0.382384 
0.365995 
0.390392 



;-2.o 

0.00 
0.10 
0.20 
0.30 
0.40 
0.50 



-0.3833072 
-0.3700764 
-0.3426137 
-0.3159575 
-0.3020421 
-0.2814511 



0.086432 
0.083773 
0.078217 
0.072782 
0.070014 
0.065737 



-0.305533 
-0.303269 
-0.297343 
-0.290017 
-0.284090 
-0.281826 



0.070364 
0.069908 
0.068718 
0.067265 
0.066138 
0.065819 



(iii) The second-order perturbation correction within the 
LF method is symmetric around ka ~ ti/2 as evident 
from Eq. (13). It is also minimum for ka = 7r/2. (iv) 
For small and intermediate values of g the second-order 
LF corrections are quite high. As a result E^ up to the 
second-order correction within the LF may not be an 
accurate estimate of Ek- Nevertheless we have cited the 
LF values for comparison with the MLF values, (v) In a 
range 1.9 < g < 2.2 (for t = oj = 1) the LF second-order 
corrections are slightly less than the corresponding MLF 
corrections. For higher values of g the MLF phonon basis 
reduces to the LF phonon basis and the results by these 
two methods become equivalent, (vi) With increasing 
g the second-order LF correction decreases considerably 

and becomes more and more /c-independent. Ef. /Ek 
within LF is around 0.027 for g — 2.5, for the whole 
region of k. 

In Fig. 3a we plot Ek vs. ka/ir for g = 0.5 and 1.0 
for the MLF approach and in Fig. 3b a comparison with 
LF result and the numerically exact results ||l2|, ^ Q 
is shown for g = 1. The MLF result is close to the nu- 
merically exact results only for low k but deviates con- 
siderably from the exact curve for higher values of k. 
In this region the second-order corrections are increas- 
ingly higher implying that higher-order corrections are 
to be included to obtain accurate results. The LF re- 
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suit is much less satisfactory than the MLF result for 
(7=1. The MLF perturbation method predicts a flat 
band (FIG. 3(b)) in the region O.Ttt < ka < tt which is 
consistent qualitatively with the exact result for 
g ~ 1.0. However, the values of Ek predicted by the MLF 
in this region is found to be not consistent with the exact 
results. 

For large g the LF second-order corrections are small 
in the entire region of fc, hence the dispersion given by the 
LF up to the second-order should be satisfactory. How- 
ever, no numerical result for Ek for large g is available to 
us for a comparison. 
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FIG. 3: Plot of Ek (in unit of t) vs ka/iv ior t ^ uj = 1.0: 
(a) Within the MLF method for g = 0.5 and 1.0. (b) Within 
the MLF and LF methods for g — 1.0; the solid line shows 
the exact results from Ref.[14]. 



In summary, in this paper we have proposed a sim- 
ple MLF phonon basis containing a single variational 
parameter but capable of describing lattice distortions 
at distant sites from the charge carrier. For low and 
intermediate values of g the second-order perturbation 
correction to the energy within this approach is much 
smaller than that corresponding to the LF approach. It 
appears that the ZPA with this MLF basis gives much 
more accurate results than that of ZPA with the LF ba- 
sis. One representative calculation of effective mass is 
presented which corroborates this feature. The ground- 
state {k = 0) energy, calculated up to the second-order 
within the MLF approach agrees well with exact results 
for the entire range of coupling strength. For g < 1.0, 
the second-order correction within the MLF is small for 
low k values and in this region the values of Ek (band 
dispersion) predicted by the MLF is consistent with the 
numerically exact result. For large g the LF perturbation 
corrections are small and the band dispersion predicted 
by the LF method is expected to be satisfactory. Com- 
putation of higher-order corrections would shed light on 
this issue. 

Electronic address for correspondence: 
moon@cmp.saha.ernet.in 
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